Карань Анна
студентка факультета биоинженерии и бионформатики

Сборка генома de novo

Сначала нужно с помощью кода доступа SRR4240356 (в моем случае) скачать файл проекта секвенирования бактерии Buchnera aphidicola в формате fastq. Эта бактерия относится к Протеобактериям и является эндосимбионтом тлей.

Этап 1

Подготовка чтений программой Trimmomatic.

java -jar /usr/share/java/trimmomatic.jar SE -phred33 SRR4240356.fastq trimm.fastq MINLEN:30 TRAILING:3 ILLUMINACLIP:adapters.fasta:2:7:7
SRR4240356.fastq - входной файл
trimm.fastq - выходной файл
MINLEN:30 - удаляет прочтения короче 30
TRAILING:3 - удаляют нуклеотиды ниже качества равного 3-м с конца прочтения
adapters.fasta - файл со всеми адаптерами для Illumina.
ILLUMINACLIP:adapters.fasta:2:7:7 - вырезает адаптеры, со значениями: 2 - отдельные несовпадения, 7 - порог для палиндромной шпильки, 7 - порог для простой шпильки.

Изначально было 7511529 ридов и файл размером 757 М, осталось 7186439 (95,67%), и размер файла - 724 М, т.е. удалено 325090 (4,33%).

Этап 2

Запуск программы velveth для подготовки k-меров.

velveth k_mer 29 -short -fastq trimm.fastq
k_mer - название папки, куда будут записаны выходные файлы
29 - длина k-меров
-short - чтения короткие и не парные
-fastq - входные файл такого формата
trimm.fastq - входной файл с очищенными чтениями.

Этап 3

Запуск программы velvetg для сборки на основе k-меров.

velvetg k_mer

И velveth, и velvetg основаны на Брюйн графах.

Таблица 1. Итоги работы программы velvetg
Финальное число узлов графа (число контигов)N50 (bp)Общая длина последовательность генома (bp) Длины 3-х самых длинных контигов (bp) Покрытия самых длинных контигов
97473133628164(ID) длина:
(8) 115468,
(1) 106076,
(9) 75082.
52,667631, 46,357791, 54,973256.

Содержимое Таблицы 1 получено при анализе файла stats.txt.
Значение типичного покрытия не так очевидно:
- 539,3478 - среднее покрытие по всем контигам
- 40,02901 - среднее контигов длиной больше 1
- 8,876263 - медиана покрытий всех контигов
- 7,186441 - медиана контигов длиной больше 1
Поэтому логичнее всего принять, что типичным покрытием является значение около 8-ми. Поэтому контигов с аномально большим покртием очень много, однако, не будем анализировать из них контиги очень маленькой длины. Вот например, контиг с покрытием 671,755906, длиной - 127, ID - 33 (контиг с самым большим покрытием длины больше 1 *под длинной один подразумевается 30 пн). На втором месте контиг с покрытием 648,875, длиной - 24, ID - 52. На третьем: покрытие - 638,8717, длиной - 413, ID - 16.
Однако, и самым низким покрытием обладают контиги длиной 1 (и не только, внизу списка по значению покрытия также короткие контиги, как и в верху). Первый более менее длиный контиг с наименьшим покрытием: покрытие - 1,5, длина - 16, ID - 415. Второй: покрытие - 1,5, длина - 20, ID - 453. Третий: покрытие - 1,578947, длина - 19, ID - 238. А наименьшее покрытие контига длины больше 1000 равно 43,61094 (ID - 29, длина - 1627), а больше 10000 - покрытие равное 46,35779 (ID - 1, длина - 106076).
Поэтому если отфилтровать контиги только длиной больше 1000 (таких 28), то максимальное покрытие будет равно 55,21522, а минимальное - 43,61094. Из этого следует, что аномальное покрытие есть только у коротких контигов! Это и логично, так как вероятность перекрытия короткого участка (или наоборот его пропуска) выше.

Этап 4

Анализ длинных контигов и контигов с аномально большим покрытием (Таблица 2).

Таблица 2. Результаты сравнения контигов с хромосомой Buchnera aphidicola с помощью megablast
ID контигаЕго особенностьКоординаты участка хромосомы, соответствующие контигу Идентичность выравнивания Число гэпов в выранивании
8Самый длинный контиг1)528794 - 550219
2)550361 - 555905
3)500370 - 508806
4)510438 - 516539
5)573092 - 582686
6)523105 - 528679
7)481997 - 488106
8)563837 - 567489
9)517766 - 521500
10)496111 - 500325
11)571558 - 573007
12)584329 - 587055
13)561741 - 563423
14)493487 - 494864
15)478095 - 480660
16)480874 - 481546
17)495033 - 495148
1)17695/21721(81%)
2)4575/5658(81%)
3)6516/8617(76%)
4)4897/6234(79%)
5)7213/9822(73%)
6)4369/5685(77%)
7)4621/6238(74%)
8)2911/3735(78%)
9)2920/3782(77%)
10)3256/4324(75%)
11)1228/1453(85%)
12)2100/2777(76%)
13)1333/1689(79%)
14)1109/1384(80%)
15)1952/2661(73%)
16)564/686(82%)
17)108/120(90%)
1)545/21721(2%)
2)133.5658(2%)
3)351/8617(4%)
4)187/6234(2%)
5)461/9822(4%)
6)207/5685(3%)
7)308/6238(4%)
8)136/3735(3%)
9)99/3782(2%)
10)154/4324(3%)
11)6/1453(0%)
12)108/2777(3%)
13)28/1689(1%)
14)13/1384(0%)
15)134/2661(5%)
16)20/686(2%)
17)5/120(4%)
12-ой по длине1)266073 - 275551
2)295935 - 303252
3)275566 - 283706
4)333222 - 339010
5)307878 - 312179
6)260224 - 263784
7)288181 - 291560
8)248967 - 252161
9)312679 - 315982
10)343228 - 346547
11)318826 - 323043
12)253244 - 257546
13)327227 - 330003
14)294227 - 295755
15)324746 - 326950
16)348233 - 349674
17)341781 - 343052
18)285200 - 286535
19)330333 - 331006
20)283963 - 285070
1)7609/9661(79%)
2)5696/7429(77%)
3)6376/8396(76%)
4)4481/5896(76%)
5)3358/4367(77%)
6)2794/3622(77%)
7)2653/3422(78%)
8)2527/3246(78%)
9)2581/3351(77%)
10)2589/3389(76%)
11)3179/4303(74%)
12)3229/4399(73%)
13)2149/2828(76%)
14)1242/1535(81%)
15)1682/2239(75%)
16)1137/1450(78%)
17)1008/1295(78%)
18)1027/1349(76%)
19)558/675(83%)
20)865/1132(76%)
1)363/9661(3%)
2)186/7429(2%)
3)421/8396(5%)
4)185/586(3%)
5)120/4367(2%)
6)111/3622(3%)
7)98/3422(2%)
8)94/3246(2%)
9)89/3351(2%)
10)118/3389(3%)
11)174/4303(4%)
12)192/4399(4%)
13)109/2828(3%)
14)14/1535(0%)
15)66/2239(2%)
16)15/1450(1%)
17)45/1295(3%)
18)27/1349(2%)
19)2/675(0%)
20)46/1132(4%)
93-ий по длине1)35124 - 44693
2)2004 - 11103
3)613658 - 620926
4)599832 - 604795
5)621055 - 627104
6)23067 - 28363
7)17962 - 20171
8)14727 - 17919
9)30028 - 32745
10)20358 - 22183
11)611633 - 613671
12)13994 - 14465
13)611229 - 611524
1)7979/9530(83%)
2)7229/9221(78%)
3)5850/7385(79%)
4)3944/5047(78%)
5)4678/6170(76%)
6)4156/5432(77%)
7)1896/2220(85%)
8)2450/3225(76%)
9)2142/2763(78%)
10)1509/1851(82%)
11)1624/2084(18%)
12)392/478(82%)
13)236/297(79%)
1)125/9630(1%)
2)252/9221(2%)
3)190/7385(2%)
4)172/5047(3%)
5)240/6170(3%)
6)217/5432(3%)
7)30/2220(1%)
8)86/3225(2%)
9)87/2763(3%)
10)51/1851(2%)
11)63/2084(3%)
12)9/478(1%)
13)2/297(0%)
SomeАномально большое покрытиеСлишком короткое выранивание, поэтому Blast его не показывает, но контигов подходящей длины и с аномально большим покрытием у меня нет.

Теперь приведем некоторые комментарии и дополнения к Таблице 2.

Рис.1. Выравнивания 8-го контига по хромосоме

Рис.2. Выравнивания 1-го контига по хромосоме

Рис.3. Выравнивания 9-го контига по хромосоме

Этап 5*

Здесь этапы 2-4 нужно проделать, поставив длину слова 25 вместо 29.

Таблица 3. Итоги работы программы velvetg при длине k-меров 25
Финальное число узлов графа (число контигов)N50 (bp)Общая длина последовательность генома (bp) Длины 3-х самых длинных контигов (bp) Покрытия самых длинных контигов
46474365628164(ID) длина:
(69) 14084,
(11) 13372,
(27) 13359.
81,245314, 81,967469, 77,609102.

Если сравнить Таблицу 1 и Таблицу 2, видно, что при уменьшении длины k-меров увеличилось число узлов (что естественно следует из строения брюйн графа), N50 значительно снизился из-за отсутствия больших контигов, самые длинные контиги в 10 раз короче для случая с k-мерами длиной 25.

Этап 6*

Здесь этапы 2-4 нужно проделать с помощью программы SPAdes.

spades.py -s trimm.fastq -k 29 --only-assembler -o spa_ass
-s - параметр, показывающий, что риды - короткие и непарные
trimm.fastq - входной файл с очищенными ридами.
-k 29 - k-меры длины 29
-o spa_ass - папка для выходных файлов

Во-первых, SPAdes работает намного дольше, чем velvet.
Ну а если о результатах, то на первый взгляд они лучше, ак как контигов меньше (62), и самый длинный из них 223925 против 115468 для velvet, разница в 2 раза. Покрытие самого длинного нуклеотида поменьше - 45,008, но это в пределах обычных значений для длинных контигов. Следующий контиг по длине - 169483, соответствено N50 - 169483, опять же в 2 раза больше, чем для velvet. Как и для velvet, аномально большие и аномально маленькие покрытия только у очень коротких контигов.
Если обобщить, то SPAdes работает лучше.

Этап 7*

Попробовать убрать половину ридов (первую половину файла) и посмотреть насколько испортится сборка.


©Карань Анна, 2015